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As a preliminary step towards simulating binary neutron star coalescing problem, we test a 
post-Newtonian approach by constructing a single neutron star model. We expand the Tolman- 
Oppenheimer-Volkov equation of hydrostatic equilibrium by the power of c -2 , where c is the speed 
of light, and truncate at the various order. We solve the system using the polytropic equation of state 
with index T = 5/3, 2 and 3, and show how this approximation converges together with mass-radius 
relations. Next, we solve the Hamiltonian constraint equation with these density profiles as trial 
functions, and examine the differences in the final metric. We conclude the second 'post-Newtonian' 
approximation is close enough to describe general relativistic single star. The result of this report 
will be useful for further binary studies. 

PACS number(s): 04.25.Dm, 04.25. Nx, and 04.40.Dg 
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I. INTRODUCTION 

Several earth-based interferometers designed to detect 
gravitational waves have been recently constructed. De- 
tectors such as LIGO, VIRGO, GEO and TAMA are ex- 
pected to begin operating within a few years (see e.g. jy). 
In order to extract gravitational waveforms from noisy 
data and to discuss physical parameters, it is essential 
to predict waveforms in advance by both analytical and 
numerical approaches. 

Binary neutron star systems are one of the most plau- 
sible sources of gravitational waves. They emit energy 
through gravitational radiation, shrink their inspiral or- 
bits gradually, and finally merge with strong emission 
of gravitational waves. The system is described by the 
post-Newtonian (PN) approximation (see e.g. Q) in the 
last several minutes before they merge, while in the last 
phase of coalescence of stars we need to solve the Einstein 
equations which are available only through numerical in- 
tegration. 

After the pioneering numerical works by Nakamura 
and Oohara in the Newtonian gravity with radiation re- 
action correction 0, several groups started developing 
their numerical codes to solve this problem in a more 
realistic way. Such hydrodynamical simulations are cate- 
gorized as in the Newtonian scheme (with/without radi- 
ation reaction term) P-|Tl|] ; PN approximation ; and 
fully general relativistic level p3| - [l5| . However, we do not 
have a method to construct physically satisfactory initial 
data for inspiral binaries in general relativity. Most of 
the numerical tests start their simulations under assump- 
tions of certain quasi-equilibrium and conformal flatness 
of spacetime, with a particular choice of vorticity of fluid 
(e.g, |ll| and references therein). 

One way to prepare initial data might be patching the 
PN scheme to the general relativistic one [|l7| . In this re- 
port, we construct a simple model and examine how this 
effort is justified. We solve the Tolman-Oppenheimer- 



Volkov (TOV) equation of hydrostatic equilibrium of a 
single neutron star, which is truncated at the various PN 
levels. We compare the mass and radius of a star as a 
function of central density using the polytropic equation 
of state. We also solve the Hamiltonian constraint equa- 
tion of the Einstein equations by substituting these den- 
sity profiles as trial functions, and discuss the differences 
in the metric. 

This study is an extended one from the earlier works 
|l8|^l| at the first PN approximation. We intend to make 
a bridge between the Newtonian and general relativistic 
solutions of a neutron star model, both of which are first 
shown numerically by Tooper [ pl[ . 

In the actual calculations, we used the geometrical 
units of c = G = Mq — 1, where c, G, M Q are the speed 
of light, Newton's gravitational constant and the solar 
mass, respectively. However, c and G will appear in the 
text where they help understanding. 



II. TRUNCATED TOV NEUTRON STARS 

In general relativity (GR), we have the TOV equation 
for solving a hydrostatic equilibrium star in the spheri- 
cally symmetric spacetime. We start from the metric 
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where e 2A M = (1- 2G e ? 2 " (r) )~ 1 . Then the TOV equations 
are written as 
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together with the specified equation of state, for which 
we use the polytropic equation of state 



p = Kp = Kp 



l+l/r. 



(2.5) 



where p, p are the pressure and energy density, respec- 
tively, and pt is the total mass density, 
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Obviously, the set of equations recover the Newtonian 
limit for c 2 — > oo. 

The idea of this rep ort is to expand the product of 
the parentheses in ( |2.3D and (fE3) and truncate them at 
the order of l/c 2 \ The i-th truncation, then, gives the 
so-called i-th PN approximation. (The case of % = 1 is 
brief ly mentioned in p3].) That is, we write (2.3) and 
(ETh schematically 
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If we use these equations with terms in the RHS of up 
to two products of A, B, C (such as AB or A 2 ), then we 
say the system is in the second PN approximation. 

We apply V = 5/3,2 and 3 for the equation of state 
(n = 1.5, 1 and 0.5 in the polytropic index, respectively) 
and compare the solutions of Newtonian, GR and up to 
third PN approximation. 

The radius of the star, R, is measured at the point, 
rv, where density p t drops low enough [O(10~ 10 ) in the 
geometrical units] , and given by the proper length, 
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with appropriate truncation in the integrand. We express 
the mass of the star, M, by M = m(r+). 

We use 5th order Runge-Kutta method (Fchlberg 
method) to integrate the equations. In order to check 
that this approach is right, we also worked the TOV equa- 
tions in the harmonic gauge and confirmed that we get 
the identical physical quantities in the results. 

In Fig|l], we show the total mass M as the function of 
the central density p c for the different Ts and PN lev- 
els. Mass is in the unit of M Q and central density is in 
[g/cm 3 ], and both are rescalable with the constant K in 
the equation of state. Here we use K in the calculations 
as: K 5/3 = 4.35 (for Y = 5/3), K 2 = 10 2 (for T = 2), and 
K3 = 10 5 (for r = 3) in the geometrical unit, where ^"5/3 
is the number for the pure neutron equation of state M . 
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FIG. 1. Total mass as the function of its central density for 
truncated neutron star model. Figures (a), (b) and (c) are for 
different equation of state with T = 5/3, 2 and 3, respectively. 
Mass is in the unit of solar mass and central density is in 
[g/cm 3 ]. The gray solid line is of Newtonian solutions, the 
solid line is of general relativistic solution. The dotted line, 
dashed line and three-dot-line are of first, second and third 
post-Newtonian approximated solution, respectively. 
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We see clearly the convergence of this PN approxima- 
tion in all the Ts. However, if the equation of state is stiff, 
then the high density configuration differs from that of 
GR even at the higher PN approximation. 

From the first PN approximation, we see the existence 
of the maximum mass. The central density which gives 
this maximum becomes larger in the weak gravity ap- 
proximation. 

In Fig.^J, we show the mass-radius relations. In the 
Newtonian limit, the asymptotic behaviors of M near 
M = are as M cx R- 3 (for r = 5/3), M cx R° (for 
r = 2) and M cx R 5 (for r = 3). These represent softness 
(for r = 5/3) and stiffness (for T = 3) of the equation of 
state. We see that all the lines in Fig. || coincide with this 
Newtonian limit in the lower mass limit. The figure also 
shows us that the first PN solution has the same feature 
as GR. 

We also checked the causality constraint dp /dp < 1 
(see e.g. 25 ) in all of the models, and confirmed that 
the constraint is always valid. 



III. METRIC OUTPUT VIA HAMILTONIAN 
CONSTRAINT 

We next solve the Hamiltonian constraint equation in 
GR with the trial density profiles obtained above. Our 
aim is to compare the difference of the output metric and 
to examine a matching scheme of PN data to the general 
relativistic one. 

We use O'Murchadha- York's conformal approach |26|] 
to solve the Hamiltonian constraint. Defining the con- 
formal factor ip and setting 7^ = V' 4 7ijj the constraint 
becomes 



(3.1) 



where and ^R are the 3-dimensional Laplacian and 
Ricci scalar curvature, respectively, defined by 7^. Here 



we assumed Kij — Kij = 0. 



(3.1 



our trial metric 7^ as conformally flat, 
with a trial density configurations of 



We choose 
and solve 

p = pt- We use the Incomplete Cholesky conjugate gra- 
dient (ICCG) method |27]] with the Robin boundary con- 
dition %jj = 1 + C/r, where C is a constant, for solving 

s 

In FigJSI, we show the conformal factor tp at the ori- 
gin as a function of central density of trial configuration. 
The 3-metric at the center will be given by 7^ = tp^dij. 
We see that using the Newotnian configuration as input 
gives us quite different solutions from the expected ones 
of GR, while all PN trials give similar solutions with GR. 
Independently to T, we can say second PN approxima- 
tion provides closer values for the output metric to those 
of GR. 
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FIG. 2. Mass and radius relations for truncated 
star models. Mass is in the unit of solar mass and 
in [km] . The lines are the same as of Fig.hl. 
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IV. DISCUSSION 
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FIG. 3. The conformal factor tjj at the origin is displayed 
as a function of central density, of which we used a trial con- 
figuration for solving Hamiltonian constraint equation. The 
central density is in the unit of [g/cm 3 ]. Each line indicates 
the trial profile as input, using the same notation with Fig.jjj. 



In order to justify the recent post-Newtonian (PN) ap- 
proaches to the binary neutron star problem, we con- 
structed a simple model. By solving the hydrostatic 
equilibrium equation of a star at i-th PN approximation, 
we showed the convergence of this approach, the mass 
and radius relations and resultant metric output via the 
Hamiltonian constraint equation. 

We conclude that second PN approximation provides 
quite similar density profiles to those of GR, independent 
of equations of states. If we use second PN density config- 
urations as trial functions, we get closer metric solutions 
to those from GR through the Hamiltonian constraint. 
Although this study is restricted to a hydrostatic single 
star model, we think that the figures shown here are con- 
venient templates for further numerical studies. 

As shown in p7| , the discontinuous matching surface 
of PN and GR in the vacuum region will be smoothed out 
in fully relativistic evolution in a particular slicing con- 
dition. Therefore we expect that higher PN initial data 
will smoothly evolve in the fully relativistic simulations, 
although there are many unknown factors as to whether 
such an initial data is numerically satisfactory or not. 
We are now applying this approach to construct a bi- 
nary model including their velocity corrections together 
with fully general relativistic hydrodynamical evolutions. 
This effort will be reported elsewhere. 
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